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We consider a system of particles undergoing the branching and annihilating reactions A — > 
(m + 1)A and A ~\- A —* %, with m even. The particles move via long-range Levy flights, where 
the probability of moving a distance r decays as r~'^~'' . We analyze this system of branching and 
annihilating Levy flights (BALF) using field theoretic renormalization group techniques close to the 
upper critical dimension dc = a, with a < 2. These results are then compared with Monte-Carlo 
simulations in d = 1. For o close to unity in d = 1, the critical point for the transition from 
an absorbing to an active phase occurs at zero branching. However, for a bigger than about 3/2 
in d = 1, the critical branching rate moves smoothly away from zero with increasing cr, and the 
transition lies in a different universality class, inaccessible to controlled perturbative expansions. 
We measure the exponents in both universality classes and examine their behavior as a function of 

O". 
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I. INTRODUCTION 

Systems possessing a continuous nonequilibrium phase 
transition from an active into an empty, absorbing state 
have been intensively studied in the past few years. De- 
spite the wide variety of processes that have been in- 
vestigated, it has proved possible to classify the criti- 
cal properties of these transitions into a small number 
of universality classes. Although the well known case 
of directed percolation (DP) |j has turned out to 
be the most common universality class, many investi- 
gations have examined systems with quite different criti- 
cal properties. For instance, the model of branching and 
annihilating random walks with an even number of off- 
spring (BARW) defines a separate universality class |^-D . 
This reaction-diffusion system consists of random walk- 
ers able to undergo the branching and annihilating re- 
actions A ^ (m + 1)A and A + A ^ 0, with m even. 
Other models in this class (at least in d = 1) include cer- 
tain probabilistic cellular automata [Q, monomer-dimer 
models [p|-p^, nonequilibrium kinetic Ising models [ pd] |, 
and generalized DP with two absorbing states . These 
models escape from the DP universality class by possess- 
ing an extra conservation law or symmetry. The BARW 
model respects an additional "parity" conservation of the 
total number of particles modulo 2. On the other hand, 
branching and annihilating random walks with an odd 
number of offspring possess no such "parity" conserva- 
tion, and hence belong to the DP universality class [^|. 
For the other models mentioned above |0-|l2[, the DP 
class is escaped via an underlying symmetry between the 
absorbing states. 

Both the DP and BARW classes do, however, share 
one important feature: the dynamical processes involved 
are short-ranged. One would expect that the addition of 
long-ranged processes would significantly alter the prop- 
erties of the active/absorbing transitions. Recently this 
expectation was confirmed by investigations of Levy DP 



(LDPV This modification, originally proposed by MoUi- 
son 1 13 1 in the context of epidemic spreading, is a gener- 
alization of DP where the distribution of spreading dis- 
tances r is given by 



P{r) - 



(a > 0), 



(1) 



where d is the spatial dimension of the system, and a is 
a free parameter (the Levy index) that controls the char- 
acteristic shape of the distribution. This distribution is 
asymptotically (as r —>■ oo) equal to a Levy distribution, 
and we will loosely refer to it as such. It was first sug- 
gested that the critical exponents describing the LDP 
transition should vary continuously with a p^ . This 
expectation was backed up by field theoretic renormal- 
ization group calculations in Ref. p5[ |, and confirmed 
numerically in Rcfs. [^,0. Note that other numerical 
work (l^Jl^ introduced an upper cut off for the flight dis- 
tance r. This resulted in effective short-range behavior, 
meaning that the LDP regime was not properly accessed. 
The results of Ref. also appear to be adversely af- 
fected by strong finite size effects. 

The purpose of the present paper is to further investi- 
gate the impact of Levy flights in models with nonequilib- 
rium phase transitions. We will analyze in detail a model 
of branching and annihilating Levy flights with an even 
number of offspring (BALF), a straightforward general- 
ization of the BARW model, where the random walkers 
are replaced by particles performing Levy flights. The 
BALF model possesses an upper critical dimension dc 
which varies continuously with the Levy index a. For 
d < dc, the model contains two new universality classes 
resulting from the long-range nature of the Levy flights. 
The exponents in both of these classes also vary continu- 
ously with a. We will investigate these new universality 
classes using field theoretic methods, some exact results 
for the pure annihilation model (where the branching pa- 
rameter is set equal to zero), and Monte-Carlo simula- 
tions in d = 1 . 
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A further attractive feature of the BALF model is that 
it casts some additional light on the properties of the 
ordinary short-ranged BARW model. We will see that 
changing the Levy index from cr = 1 to cr 2 for the 
BALF model in fixed dimension d = 1 is in many ways 
similar to changing the physical dimension from d = 2 to 
d = 1 in the short-ranged BARW model. Although this 
correspondence is certainly not rigorous, we can never- 
theless use simulations of the BALF model in the physical 
dimension d = 1 to better understand properties of the 
BARW model which lie in the inaccessible dimensions 
between d — 1 and d = 2. This will allow us to probe 
numerically some important features of the BARW field 
theory developed by Cardy and Tauber in Ref. ||] . 

We now give a brief summary of the layout of the pa- 
per. In the next section we briefly review the relevant 
properties of the short-ranged BARW model. In Sec- 
tion III we then introduce the BALF model and present 
its mean field behavior. After these preliminaries we then 
present the field theoretic action for BALF, which we 
analyze using diagrammatic and renormalization group 
methods. These results are then compared with Monte- 
Carlo simulations in Section IV. Finally, our conclusions 
appear in Section V. 



II. BRANCHING AND ANNIHILATING 
RANDOM WALKS 

The BARW model is defined by the following reaction 
processes: 



A 



{m+l)A, 



rate firn, 
rate A, 



(2) 
(3) 



where the identical particles A otherwise perform simple 
random walks with diffusion constant D. As the reac- 
tion rate parameters are varied, one finds a continuous 
phase transition from a region controlled by the pure an- 
nihilation process to an active region characterized by a 
nonzero particle density in the steady-state. The growth 
of BARW clusters close to the critical point can be sum- 
marized by a set of independent exponents. A natural 
choice is to consider and ;/|| , which describe the diver- 
gence of the correlation lengths in space, ^± ^ |A|~'^-l, 
and time, ^|| ~ lAI""^!!, close to criticality. Here the pa- 
rameter A describes the deviation from the critical point 
at the active/absorbing transition. We also need the or- 
der parameter exponent /3, which can be defined in two 
a priori different ways: it is either governed by the prob- 
ability that a cluster grown from a finite seed never dies. 



P{t ^ oo, A) - A^- 



A > 0, 



(4) 



or by the coarse-grained density of active sites in the 
steady state. 



n(A) ~ A^-* 



These exponents can be simply calculated in mean field 
theory, valid for d > d^ — 2. The appropriate mean 
field rate equation for the coarse-grained density n(x, t) 
is given by 



mfimn — 2\r? . 



(6) 



For /i„i = no branching is present, and we are reduced 
to the well known annihilation reaction A + A 0, which 
asymptotically exhibits a power law mean field density 
decay n ~ t^^. However, for nonzero /i™, we have the ho- 
mogeneous steady-state solution = m/i,„/2A. Hence 
the critical value of clearly lies at zero, and we iden- 
tify A = m^m- The density thus behaves as ris oc A, 
and we immediately see that = 1. The alternative 

order parameter exponent f3^^ can also be simply cal- 
culated: for d > dc — 2, the survival probability of a 
particle cluster will be finite for any value of the branch- 
ing rate, implying that = 0. This result follows 
from the non-recurrence of random walks in d > 2. The 
correlation length exponents can also be simply derived 
from Eq. (^), yielding i^^^ = 1/2, i^h^^ — 1. Hence the 



dynamic exponent z, defined by ^ ^1/^^ given by 



,MF /,,MF 



Below the upper critical dimension the above mean 
field analysis breaks down due to the presence of fluctu- 
ations. Recently, methods have been developed to sys- 
tematically include these fluctuation effects. Firstly, the 
appropriate master equation, which provides a complete 
description of the microscopic dynamics of the system, is 
transformed into a second-quantized Hamiltonian. This 
representation is then mapped onto a coarse-grained field 
theoretic action p0|-|2^. From this point the standard 
tools of renormalized perturbation expansions can be em- 
ployed, and the effects of fluctuations systematically com- 
puted. For the case of BARW, this analysis was per- 
formed in Ref. Q. In the following we summarize the 
main results of that analysis. The field theoretic action 
for BARW, written in terms of the response field t) 
and the "density" field 'ijj{-x.,t), is given by P| 



5'o[V','0;t] 



d^x 



dt 



V'(x,t)[a- W2]^(x, t) 
A[l-^(x,i)2]V(x,t)2 (7) 



A > 0. 



(5) 



+/x™[l-^(x, i)™]?/.(x,i)V'(x,t) 



-V'(x, r) - no^(x,0) 



Here the terms on the first line of (|^) represent diffu- 
sion of the particles (with continuum diffusion constant 
D). The second line describes the annihilation reaction 
(with continuum rate A), while the terms on the third 
line represent the branching process (with continuum 
rate Hm)- The final two terms represent, respectively, 
a contribution due to the projection state (see Ref. pof). 
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and the initial condition (an uncorrelated Poisson distri- 
bution with mean ng). In the following we will restrict 
ourselves to the case of even m, since it is known that 
the odd m case belongs to the DP universality class . 

The action given in (Q) is a bare action. In order to 
properly include fluctuation effects one must be careful to 
include processes generated by a combination of branch- 
ing and annihilation. In other words in addition to the 
process A ^ {m + 1)A, the reactions A ^ {m— I) A, . . . , 
A 3A need to be included. These considerations lead 
to the full action 



n/2 



dt 



i/;(x,t)[at-Dv2]V'(x,i) 

^M2i[l-V'(x,t)'']^(x,i)V'(x,i) (8) 



1=1 



-A[l - V3(x, i)2]i^(x, tf ~ V^x, r) - no7/^(x, 0) 



Notice also that (for even m) the action (^ is invariant 
under the "parity" transformation 



V'(x,t) -V'(x,t), 



V;(x,<) 



-V^(x, t). 



(9) 



This symmetry corresponds physically to particle con- 
servation modulo 2. The presence of this extra symme- 
try now takes the system away from the DP universality 
class, and into a new class: that of branching and anni- 
hilating random walks with an even number of offspring. 

Simple power counting on the action in Eq. (H) reveals 
that the upper critical dimension is dc = 2. Close to dc, 
the renormalization of the above action is quite straight- 
forward (here we again quote the results from Ref. [^). 
At the annihilation fixed point the RG eigenvalue of the 
branching parameter can easily be computed. To one 
loop order one finds y^^^ = 2 — m{m + l)e/2 + O(e^), 
where e = 2 — d. Hence we see that the lowest branching 
process is actually the most relevant. Therefore, close to 
2 dimensions where the branching remains relevant, we 
expect to find an active state for all nonzero values of 
the branching (in agreement with the mean field theory 
presented above). Furthermore, in this regime, we can 
exploit the fact that the critical point, which remains 
at zero branching, is described by the pure annihilation 
theory. Matching the exactly known density decay pl| ] 
and survival probability exponents in the annihilation 
theory with their counterparts in the critical BARW the- 
ory yields the exact exponent relations /3dons = dz/||/2, 
Psccd = J^||(2 — d)/2, and z = 2 ||]. To the best of our 
knowledge, the result for /?socd, although simple to derive, 
has not previously been given in the literature. 

Inspection of the one loop result for the most rele- 
vant RG eigenvalue shows that it eventually becomes 
negative. This occurs at a second critical dimension d^, 
where d^ = 4/3 to one loop order. For d < d^ we expect 
a major change in the behavior of the system, since the 
branching process will no longer be relevant at the anni- 
hilation fixed point. The critical transition point is then 



shifted with the active state only being present for values 
of the branching greater than some positive critical value. 
For branching parameter values smaller than this value, 
the branching will be asymptotically irrelevant. This re- 
gion of parameter space will thus be controlled by the 
annihilation fixed point of the A + A^% process, where 
the density decays away as a power law. Hence this re- 
gion of parameter space should be considered a critical 
inactive (or absorbing) phase. The presence of a second 
critical dimension d^ « 4/3 immediately rules out any 
possibility of using perturbative e expansions to access 
the non-trivial active/ absorbing transition expected in 
the physical dimension d = 1. Instead cruder techniques 
(such as the loop expansion in fixed dimension) must be 
employed |^ . We will not discuss this part of the analysis 
of Ref. 1^ in much detail. However we do wish to point 
out that the truncated loop expansion at one loop does 
predict a jump in the critical point at around d^ « 4/3, 
from zero branching to some finite value. We will have 
more to say about his observation in Section IV, after we 
have presented our analytical and numerical study of the 
BALF model. 

It is also possible to analyze the BARW model in d = 1 
using exact methods. In Ref. |Q it was demonstrated 
that, at the annihilation fixed point, the one loop RG 
eigenvalue y^J_^ is actually exact in d — 1. However, 
the reason for the cancellation of the contributions from 
higher loop orders in the field theory remains unclear. 
Other work |23|j2^ ], using quantum spin Hamiltonians, 
has indicated that the exponents /?dons and /?scod are ex- 
actly equal at the active/ absorbing transition in d = 1. 
This conclusion is also supported by numerical simula- 
tions l25tl. 



III. BRANCHING AND ANNIHILATING LEVY 
FLIGHTS 

We now turn to the main object of this paper: a sys- 
tematic investigation of the BALF model. To begin with, 
we consider the model at the mean field level. The ap- 
propriate mean field equation is given by 



dtn = (DatV^ + DA^")n + mii^n - 2\t? 



(10) 



where and Da are the rates for normal and anoma- 
lous (Levy) diffusion, respectively. The anomalous dif- 
fusion operator describes moves over long distances 
and is defined by its action in momentum space 



(11) 



where k — |k|. The standard diffusion term Dm^^ takes 
into account the short-range component of the Levy dis- 
tribution. A more detailed derivation and justification 
for the Levy term can be found in Ref. 0- The mean 
field exponents can now easily be extracted. The criti- 
cal point remains at zero branching, and, for cr < 2, we 
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jMF 
seed 



0, i^f ^ 



1, and = l/a. 



identify = h A 

Note that, for cr > 2, these exponents cross over smoothly 
to the ordinary mean field BARW exponents. Even at the 
mean field level, we can see that the exponent v± varies 
continuously with the Levy index cr. 

The above mean field description will only be quan- 
titatively valid above the upper critical dimension. For 
d < dc, we must again take fluctuation effects into ac- 
count. This can be done using the same methods as were 
used for the short-ranged BARW model |^ . We empha- 
size that the inclusion of the long-ranged Levy processes 
does not introduce any particular difficulties for the field 
theory mapping (see Ref. jl^ for further details). Spe- 
cializing immediately to the case with m — 2, and defin- 
ing /i = /^2 7 we find that the field theoretic action is given 
by 



d'^x 



dt 



i/;(x, t)[dt - DnV^ - i?^V'^]^(x, t) 

-A[l-^(x,i)2]V;(x,t)2 (12) 
+Ai[l-V3(x,i)2]i^(x,t)V(x,t) 

-V'(x, r) - no'(/'(x, 0) . 



This action describes both normal and anomalous diffu- 
sion. The naive scaling dimensions of the fields are 



[V;(x,i)] =^0, [V'(x,i)]=K^ 



(13) 



With [x] — K ^ and [<] = k we see that the naive 
scaling dimensions of the couplings are 



(14) 



Hence, power counting reveals that the upper critical di- 
mension, at which the fluctuations become important, is 
dc = cr, for cr < 2. 

We have calculated the renormalization group flow 
functions and eigenvalues, so as to determine the long 
distance and late time behavior of this field theory. The 
one loop contribution to the renormalized annihilation 
vertex is given by the diagram in Fig. |^a. For the case 
where n = 0, the propagator is (s + D^k'^ + DAk'^)^^ in 
(k, s) space (s is the Laplace transformed time variable), 
or ^-{DAk'+Dnk )t £^ space. It turns out to be easi- 

est to calculate an extended-time vertex function in (k, t) 
space, and then determine the renormalized coupling in 
(k, s) space by performing a Laplace transform and eval- 
uating at the normalization point (k, s) = (0, 2Z)^k°'). 



(a) 





(b) 



FIG. 1. One loop contribution to (a) the renormalized an- 
nihilation vertex and (b) the renormalized branching vertex. 

The first step is to drop the normal diffusion term, as 
it is less relevant for cr < 2. The dimensionless renormal- 
ized annihilation coupling is defined by 



with e = dc — d = (7 — d, and 

r{d/a) r(2 - d/a) 



Cd 



T{d/2) 2d 



-l_-d/2 



The one loop renormalization factor Zx is then 

Zx = i . 

Da e 

Hence the (3 function is given by 



(15) 



(16) 



(17) 



i{d- 



(18) 
a -d. The 



with fixed points at ^ = and I — i* = 
result (18) is actually exact to all orders in perturbation 
theory^. For d > cr, the Gaussian fixed point at ^ = 
is stable, while for d < a, the non-trivial 0(e) fixed point 
a.t £ = £* is stable. 

To investigate the relevance of the branching process, 
we now calculate the one loop RG eigenvalue for the 
branching process at the annihilation fixed point. Defin- 
ing the dimensionless renormalized branching rate as 



Z^fiK " /Da, 



(19) 



then, from the diagram in Fig. we can compute the 
one loop renormalization factor 



Z„ 



1 



, A CdK 
'Da e 



Hence the C, function is 
d 



OK Id 



3£ + Oi£^ 



(20) 



(21) 



Thus the one loop RG eigenvalue for the branching pro- 
cess at the annihilation fixed point is 



y^, = -C^,i£*) = (J-Se = Sd-2a. 



(22) 
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Consequently, according to the one loop theory, the 
branching process is relevant at the annihilation fixed 
point for 

a < a'^{d) = id/2, (23) 

or, in d = 1, for a < a'^{d = 1) = 3/2. Hence, as in 
the mean field case, we expect an active phase for all 
nonzero values of the branching rate fj., for sufficiently 
small a (see also the phase diagram in section IV). In 
this regime, we can again exploit the fact that criticality 
lies at zero branching, and hence that the critical behav- 
ior of the BALF model coincides with that for the simple 
Levy annihilation model A + A 9 [jl^ . For the Levy 
annihilation model, several exact results can be derived: 
the density decays as t~'^^'^ (for d < a < 2) the sur- 
vival probability decays as t'^/<^~^ (also for d < a < 2); 
and the dynamic exponent is just z ^ a (for a < 2). 
The second of these results follows in a simple way from 
the analysis of Ref. |^6|, but is nevertheless, to the best 
of our knowledge, a new result. On the other hand in 
the BALF model it is straightforward to show that, at 
criticality, the density should decay as i-Z^dens/wy ^ 
the survival probability as t~''=<=<=d/'^ii ^ and where again 
we have z = cr at the Levy annihilation fixed point p7[ . 
Matching these results to the Levy annihilation case, we 
have /3dons = dv\\la and /3sccd = v\\{cr - d)/(j. 

Hence, in the regime a < cr'^{d), there is just one in- 
dependent exponent which must be calculated perturba- 
tively. Following a similar analysis to that in Ref. |Q , this 
exponent can be taken to be i'± , which in terms of the 
RG eigenvalue for the branching is given by i'± = 
Ideally, at the annihilation fixed point in d = 1, one 
would like to be able to calculate the RG eigenvalue 
exactly, as was done in the short-ranged BARW model 
1^. Remarkably in that case it was found that the one 
loop result was exact. Unfortunately, generalizing the 
methods of Ref. |^ to the Levy case does not seem to 
be straightforward. However there is some numerical ev- 
idence (to be presented in the next section) to suggest 
that the one loop result for in the Levy case is again 
exact. 

We also note that, in the above regime, if a = 2 < a'^ 
then, to one loop order, we expect a smooth crossover to 
the short-ranged BARW model. Hence, at least at the 
one loop level, the model in this regime is more straight- 
forward than the LDP case^where there are additional 
complications (see Refs. ||l5| , [l7|| for more details). 

We now discuss the case where, for a > cr'^{d), the 
branching becomes irrelevant at the Levy annihilation 
fixed point. For this regime to be present at all, then 
from Eq. (p3|), we require d < 4/3 to one loop order. In 
this regime we expect the critical branching rate fJ-dcr, d) 
to become nonzero. For < < /Xc((t, cZ), the branch- 
ing will be asymptotically irrelevant, and this phase will 
again be governed by the exponents of the pure Levy an- 
nihilation universality class p7[ |. At fi = i^do', d) we then 
expect a non-trivial transition to an active phase. As was 



the case for the short-ranged BARW model we expect 
this transition to be inaccessible to controlled perturba- 
tive expansions, and in a different universality class to 
that discussed above. This follows from the fact that this 
transition only appears below a second critical dimension 
d < (ic(cr) = 2(t/3 to one loop order. Hence, as was the 
case for the short-ranged BARW model, e = dc — d ex- 
pansions down from the upper critical dimension dc = u 
will not be able to access this transition. 

We note that a "precursor" of the critical inactive 
phase present for a > cr'dd) is already evident in the 
a < cr'c{d) regime. Using the above analysis we see that 
/3dons = d/y^ = d/{3d — 2a) for a < (j'^d) (to one loop). 
Hence, to this order, as cr ^ <^'c.{d) = id/ 2 from below, 
/3dons diverges. This implies that for a fixed small value 
of the branching, the active phase has a decreasing den- 
sity as a function of a, as a is increased towards cr'dd). 
Finally at cr = 3d/2 (to one loop order) an entire critical 
inactive phase opens up. 

We can now see the similarities between properties of 
the short-ranged BARW model as dimension is lowered 
from d = 2 to c? = 1, and the BALF model 'md — 1, as the 
Levy index is raised from a = 1 to a = 2. In particular, 
to one loop order, the region 1 < cr < cr'dd = 1) = 3/2 
for the d = 1 BALF model contains the direct analog of 
the inaccessible universality class present in BARW for 
< « 4/3 < d < 2. 



IV. SIMULATION RESULTS 

In order to further investigate the BALF universality 
class, we have performed extensive numerical simulations 
of a lattice BALF model in d = 1. At each time step, a 
randomly chosen particle was allowed either to branch, 
with probability 1 — p, or to move, via a long-range jump, 
with probability p; p was the only parameter in the sim- 
ulations. The number of particles at each lattice site was 
restricted to zero or one: thus, when a particle moved to 
an occupied site, both particles were annihilated. At each 
branching step, a particle produced two offspring, which 
occupied the two sites to the immediate left or right of 
the original particle, with the side chosen randomly. As 
pointed out in Ref. |^^, this method of choosing occu- 
pied sites is necessary since if the newly occupied sites 
are chosen symmetrically about the original site, then 
the short-ranged BARW model turns out to be in its 
inactive state for all < p < 1. 

The distribution of hop lengths was chosen to follow 
Eq. for r > 1. This distribution is implemented by 
choosing a random number x from a uniform distribution 
on the interval [0, 1), and then calculating a new random 
variable r = (1 — x)"^/"^. It is easy to see that this pro- 
duces a sequence of numbers whose distribution follows 
Eq. (0). 

Two different initial conditions were used to calculate 
the different exponents. In the first case, the initial con- 
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dition was a "seed" of two particles at lattice sites ±1. 
See Fig. ^ for some sample runs with this initial condi- 
tion, run for 500 time steps. The long-range hops at 
small a result in a very rapid and wide-ranging dispersal 
of the particles. 



-4x10 



4x10 -4(X) 




4(X) 



200 



200 







200 



-■2(X) 

FIG. 2. Sample runs at various values of a, with time evo- 
lution running up the page. AH runs are for 500 time steps, 
and are at values of p about 10% away from the critical point 
into the active phase. The top two frames show p = 0.9, at 
a = 1.0 and a — 1.5, from left to right. The lower left frame 
shows a = 1.9, p = 0.77, and the lower right frame shows the 
ordinary short-ranged BARW model at p = 0.46. Notice the 
large change in scale between the first and last frame. 

These simulations were averaged over many runs from 
the same initial condition but for different sequences of 
random numbers. The number of runs, P{t), surviving 
to time t, the number of particles in the system averaged 
over the total number of runs, N{t), and a mean square 
spreading distance, R^{t), were all measured. The mean 
spreading distance is defined by a geometric mean in the 
Levy case (see Ref. [0 for more details). At the crit- 
ical point, these quantities should all follow power law 
behavior with 



Pit)' 
N{t) 



2/z 



(24) 
(25) 
(26) 

The critical point is determined by plotting a local ex- 
ponent against 1/t, and estimating the value of p which 
produces a straight line as 1/t 0. For the survival 
probability, the local exponent is defined by 



m 



In 



Pjt) 



In 6 



(27) 



and similarly for the other quantities. We have used & = 5 
in our data analysis. The extrapolation of the exponent 
to 1/i ^ is the estimate of its long time value. A 
sample analysis, for a = 1.6, is shown in Fig. 0. 



2/z 




FIG. 3. The effective local exponents, as in Eqs. (|24|)-(g6|), 
for a = 1.6. The exponents are plotted against 1/t to extract 
the t — + oo limit. The curves correspond to values of p, from 
top to bottom, of 0.985, 0.988, 0.99, and 0.992. 

We focus first on the perturbatively inaccessible tran- 
sition found for larger values a > <j'^{d = 1). In this 
regime we performed simulations at criticality and mea- 
sured the exponents defined in Eqs. (24)-(p6[). To reduce 
finite size effects, we implemented periodic boundary con- 
ditions and used a very large lattice. Rather than store 
the occupation numbers of each lattice site, only the po- 
sitions of the particles were stored. This meant that the 
system size was limited only by the number of integers, 
i.e. a system size of 2^'^ w 1.8 x 10^^ on the 64-bit com- 
puter used. The simulations ran for times of between 
2 X lO** and 2 x 10^ time steps per particle, and were 
averaged over at least 2 x 10® runs. 

We encountered several obstacles in accurately deter- 
mining the values of these exponents. First, the quan- 
tities measured were also expected to behave as power 
laws on the (critical) inactive side of the transition. The 
difference between the equivalent exponents on the crit- 
ical line and in the critical phase were sometimes small, 
particularly near a = (j'dd = 1), where the fixed point 
we investigated merges with the pure Levy annihilation 
fixed point. Consequently measurements close to this 
point required the longest runs. Also, corrections to scal- 
ing meant that the effective local exponents did, in fact. 
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vary with t. The impact of these corrections to scahng 
was sometimes difficult to interpret accurately. Finally, 
the various exponents measured yielded slightly different 
estimates for pc, making it difficult to determine the crit- 
ical point more accurately than done here. This may be 
due to the corrections to scaling being of different sizes 
for each of the exponents measured. 

The results of these simulations are shown in Table |. 
The fj = oo values shown are for simulations with normal 
diffusion, where the exponents measured are consistent 
with those of other simulations 27|. Fig. |^ shows the 
phase diagram as determined by the simulations. 

We now discuss some features of the numerical data in 
Table |: 

• The data presented in Table || is consistent with a 
value very close to cr = a'^{d = 1) = 3/2 for the emer- 
gence of the critical Levy annihilation phase at nonzero 
branching. This is i n g ood agreement with the one loop 
result for j/^ in Eq. ( p2|) , and provides some evidence that 
this one loop result may in fact be exact, as it was for 
the short-ranged Levy case. 

• The measured exponents changed by rather small 
amounts over the range of a studied. As discussed in 
section III, the exponents at cr = a'^{d = 1), p = 1, can 
be calculated for the pure Levy annihilation model, and, 
assuming a'^{d = 1) = 3/2, are given hy 6 = —6 = 1/3, 
z = 3/2. If the exponents are to change monotonically 
as a is varied, then they are trapped in a relatively small 
range of values between the BARW and Levy annihila- 
tion exponents. 

• The numerical evidence also strongly favors a smooth 
movement of the critical value Pc away from unity as 
cr is increased above 3/2, as shown in Fig. ^. This 
finding has consequences for the analogous short-ranged 
BARW model. In that case the analog of the point at 
cr = = 1) w 3/2 is the second critical dimension 
found at d'^ w 4/3. The uncontrolled truncated loop ex- 
pansion used to analyze this point in Ref. |^ predicted 
a discontinuous jump of the critical point as dimension 
d is lowered through d'^. The above numerical evidence 
argues against this scenario and would rather predict a 
smooth movement of the critical point. Given the uncon- 
trolled nature of the truncated loop expansion, any fail- 
ure to accurately capture the behavior close to d'^ would 
not, perhaps, be very surprising. Nevertheless, our re- 
sults have, for the first time, provided numerical evidence 
for one of the main conclusions of Ref. namely the 
presence of a second critical dimension d'^. 

• Despite considerable effort, the data reported in Ta- 
ble I are unfortunately not precise enough to answer the 
question: at what value of a do the Levy results cross 
over to those of the short-ranged BARW model? Re- 
grettably, the situation from a theoretical perspective is 
no clearer, due to the absence of any controlled field the- 
oretic methods in this regime. 

We now turn our attention to the second regime for 
d=l BALF, that for a < a'^{d = 1) w 3/2. In this case, 
it is not appropriate to perform simulations at critical- 



ity, since in that case we would only be measuring the 
exponents of the pure Levy annihilation model. Hence 
we have performed off-critical simulations in an effort 
to measure /?dcns as a function of cr. In this case, we 
used a second initial condition, a fully occupied lattice 
of size L — 10^. We then allowed the number of parti- 
cles to decay away until a steady-state was reached. The 
steady-state density depends on the deviation from the 
critical point, as described by Eq. (^, and thus /3dcns may 
be directly measured. 

/ Critical Annihilation Phases 



P 0.8 



0.7 



0.6 







Active Phase \. \i 

1 \ 








0.95 








1.5 ^ 1.7 ^ \ 



1 



1.5 



2.5 



FIG. 4. Phase diagram for the BALF model in d = 1. The 
inset shows a blowup of the region near a'c{d — 1). 



a 


Pc 


5 


e 


z 


1.525 


0.997(2) 


0.32(1) 


-0.30(1) 


1.53(2 


1.55 


0.992(4) 


0.32(2) 


-0.30(2) 


1.53(2 


1.6 


0.990(2) 


0.33(2) 


-0.30(2) 


1.56(2 


1.65 


0.974(2) 


0.32(2) 


-0.26(2) 


1.55(2 


1.7 


0.955(2) 


0.32(2) 


-0.24(2) 


1.59(2 


1.8 


0.918(2) 


0.32(2) 


-0.18(2) 


1.59(2 


1.9 


0.863(2) 


0.32(2) 


-0.14(1) 


1.63(2 


2.0 


0.804(1) 


0.305(5) 


-0.085(5) 


1.68(2 


2.5 


0.6185(2) 


0.285(5) 


-0.005(5) 


1.72(1 


oo 


0.5104(2) 


0.287(3) 


0.001(3) 


1.74(1 



TABLE I. The measured critical probabilities and expo- 
nents for various values of a. The number in brackets is an 
estimate of the error in the last figure. 
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The values of /3dcns measured in these steady-state sim- 
ulations are given in Table |l|. For ct < 1 the mean field 
result should hold, since d = 1 then lies above the upper 
critical dimension dc — o. For cr slightly bigger than 
unity, the upper critical dimension will lie just above 
d — \^ and hence one might hope to directly observe 
the e expansion results (see also Ref. for a similar 
case) . Unfortunately, the values measured in the simula- 
tions deviate by around 5 — 25% from the mean field or 
one loop e expansion exponents calculated in section III. 
We believe there are two reasons for this discrepancy. 
Firstly, for small ct, finite size effects become important, 
as the long-ranged hops allow a single particle to wrap 
all the way around the system in a short time. Secondly, 
as CT — > a'Jyd = 1) w 3/2, /3dcns becomes rather large, 
and hence large systems and long runs were necessary to 
probe the very small steady-state densities which occur 
near pc- Although we used as large a system as practica- 
ble, we were not able to entirely eliminate the discrepancy 
between theory and simulations. 

In summary, despite the difficulties encountered for 
CT < o'J^d — 1), the overall picture that emerges from 
the numerics agrees well with the theory presented in 
the last section. As we have discussed earlier, this, in 
turn, provides additional support for the analysis of the 
short-ranged BARW model presented in Ref. 



a 


/3dcns (measured) 


/Jdens (theory) 


0.7 


1.0 


1 (mean field) 


0.9 


1.1 


1 (mean field) 


1.1 


1.3 


5/4 (one loop) 


1.3 


1.8 


5/2 (one loop) 



TABLE II. The exponent /3dens determined in simulations 
of the steady-state density of a system of size 10^. 



V. CONCLUSIONS 

In this paper we have presented an analytic and nu- 
merical study of the BALF model. Using field theoretic 
techniques we have obtained a good analytic understand- 
ing of the model in the physical dimension d = \ for the 
regime ct less than about 3/2. For values of ct larger 
than this, we have had to rely solely on numerical sim- 
ulations. In both regimes the critical exponents of the 
active/absorbing transition are found to vary continu- 
ously with the Levy index ct. Numerically, we find that 
the transition between the two regimes in d = 1 occurs 
at CT = CT^(d = 1) « 3/2, in agreement with the one loop 
result from the field theory. Unfortunately our numerics 
for the small ct regime were not good enough to confirm 
the accuracy of our e expansion calculations. Neverthe- 
less this is the first time this universality class has been 
accessed, since its equivalent in the original BARW model 
lies in the inaccessible dimensions < < 2. 

Finally, we would like to emphasize that Levy flights 
are a powerful way of probing the higher dimensional be- 
havior of nonequilibrium models while performing simu- 
lations only \x\ d— 1 . The disadvantage of this approach 
is that it necessitates the use of extremely large system 
sizes if finite size effects are to be avoided. However, we 
have shown that in large regions of parameter space these 
problems can be overcome, and reasonable estimates ob- 
tained for the exponents. 
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